library(tidyverse)
library(skimr)
library(patchwork)
#library(maps)
library(janitor)
library(lubridate)
library(leaflet)
library(mapdeck)
theme_set(theme_classic())Basic_LOC01_plots
Basic LOC01 plots from Jennie’s python notebook, translated to R.
Data
Load data. Files from LOC-01 data on Teams drive.
loc01_underway <- "data/LOC-01_Underway_continuous.txt"
loc01_drifter_comp <- "data/drifters_interp.csv"
loc01_ctd <- "data/CTD_downcast_upcast.csv"Underway
uw_df <- read_csv(loc01_underway,
na = c("", "NA", "NaN"),
col_types = "cnninnnnnncni",
) |>
mutate(DateTime_UTC = as.POSIXct(DateTime_UTC,
format = "%d-%b-%Y %H:%M:%S",
tz = "UTC"),
Salinity_PSU = ifelse(Salinity_PSU < 29, NA, Salinity_PSU))
uw_dye_df <- uw_df |>
filter(DateTime_UTC > as.POSIXct("2023-09-02 12:00:00", tz = "UTC"),
DateTime_UTC < as.POSIXct("2023-09-04 01:30:00", tz = "UTC"))CTD
ctd_df <- read_csv(loc01_ctd) |>
clean_names()Need per-cast info: station, cast, in-out, start timestamp, start lat, start lon
ctd_loc <- ctd_df |>
mutate(timestamp = as.POSIXct(paste(start_date, as.character(start_time_utc)), format="%m/%d/%y %H:%M:%S")) |>
select(station, cast, in_out,
timestamp,
lat = longitude_deg, lon = longitude_deg_2) |>
distinct(station, cast, .keep_all = TRUE)Drifters
Data for drifters 3 and 4 processed in separate script. Position and time between SPOT fixes are interpolated (spline). One outlier in drifter 3 and some position data without sonde data in drifter 4 removed.
drifter_df <- read_csv(loc01_drifter_comp)Timeseries
Underway
dye_plot <- ggplot(uw_dye_df, aes(DateTime_UTC, Dye_ppb)) +
geom_line() +
scale_y_log10() +
labs(x = "Time",
y = "Dye (ppb)")
dye_plot_ctd <- ggplot(uw_dye_df, aes(DateTime_UTC, Dye_ppb)) +
geom_point() +
geom_line() +
scale_y_log10() +
labs(x = "Time",
y = "Dye (ppb)")
temp_plot <- ggplot(uw_dye_df, aes(DateTime_UTC, Temp1_C)) +
geom_line() +
labs(x = "Time",
y = "Temp (C)")
sal_plot <- ggplot(uw_dye_df, aes(DateTime_UTC, Salinity_PSU)) +
geom_line() +
labs(x = "Time",
y = "Salinity (PSU)")
ta_plot <- ggplot(uw_dye_df, aes(DateTime_UTC, Corrected_TA_umol_kg_)) +
geom_point(shape = 1) +
labs(x = "Time",
y = expression(paste(TA ~ (μmol ~ kg^-1))))
(dye_plot + temp_plot) / (sal_plot + ta_plot)
CTD Plots
#|fig-height: 800
ggplot(ctd_df, aes(dep_sm_salt_water_m, rhodfl_tc0_ppb, color = upcast_downcast)) +
geom_line() +
coord_flip() +
scale_y_log10() +
scale_x_reverse() +
facet_grid(rows = vars(in_out), cols = vars(station, cast))
Drifters
drifter_df |>
ggplot(aes(datetime, dye, color = asset)) +
geom_point() +
scale_y_log10()
drifter_df |>
ggplot(aes(datetime, dye)) +
geom_line(data = uw_dye_df, aes(DateTime_UTC, Dye_ppb), color = "lightblue") +
geom_line() +
scale_y_log10() +
facet_grid(rows = vars(asset))
Maps
Underway
# Outline map
map_data = map_data("usa")
myscaler <- function(x, from, to) {
high=1.2
low=0.1
ifelse(x<low,0,ifelse(x>high,1,(x-low)/(high-low)))
}
# Plot
ggplot() +
geom_map(data = map_data, map = map_data,
aes(long, lat, map_id = region),
fill = "white", color = "black") +
geom_point(data = uw_dye_df,
aes(x = Longitude, y = Latitude, color = log10(Dye_ppb)),
size = .4) +
scale_color_viridis_c(option = "D",
name = "Dye (log ppb)",
#trans = scales::log10_trans(),
rescaler = myscaler) +
xlim(-70.8, -70.6) +
ylim(41.01, 41.15) +
coord_quickmap() +
ggtitle("Dye") +
xlab("Longitude") +
ylab("Latitude")
ggplot() +
geom_map(data = map_data, map = map_data,
aes(long, lat, map_id = region),
fill = "white", color = "black") +
geom_point(data = uw_dye_df,
aes(x = Longitude, y = Latitude, color = Salinity_PSU),
size = .4) +
scale_color_viridis_c(option = "D") +
xlim(-70.8, -70.6) +
ylim(41.01, 41.15) +
coord_quickmap() +
ggtitle("Salinity") +
xlab("Longitude") +
ylab("Latitude")
ggplot() +
geom_map(data = map_data, map = map_data,
aes(long, lat, map_id = region),
fill = "white", color = "black") +
geom_point(data = uw_dye_df,
aes(x = Longitude, y = Latitude, color = Temp1_C),
size = .4) +
scale_color_viridis_c(option = "B",
name = "Temp (C)") +
xlim(-70.8, -70.6) +
ylim(41.01, 41.15) +
coord_quickmap() +
ggtitle("Temperature") +
xlab("Longitude") +
ylab("Latitude")
cruise_TA <- uw_df |>
filter(!is.na(Corrected_TA_umol_kg_)) |>
ggplot() +
geom_map(data = map_data, map = map_data,
aes(long, lat, map_id = region),
fill = "white", color = "black") +
geom_point(aes(x = Longitude, y = Latitude, color = Corrected_TA_umol_kg_),
size = 2) +
scale_color_viridis_c(option = "D",
name = expression(paste(TA ~ (μmol ~ kg^-1)))) +
xlim(-72.1, -70.4) +
ylim(40.5, 41.5) +
coord_quickmap() +
ggtitle("Alkalinity") +
xlab("Longitude") +
ylab("Latitude")
dye_TA <- uw_dye_df |>
filter(!is.na(Corrected_TA_umol_kg_)) |>
ggplot() +
geom_map(data = map_data, map = map_data,
aes(long, lat, map_id = region),
fill = "white", color = "black") +
geom_point(aes(x = Longitude, y = Latitude, color = Corrected_TA_umol_kg_),
size = 2, shape = 1) +
scale_color_viridis_c(option = "D",
name = expression(paste(TA ~ (μmol ~ kg^-1)))) +
xlim(-70.8, -70.6) +
ylim(41.01, 41.15) +
coord_quickmap() +
ggtitle("Alkalinity") +
xlab("Longitude") +
ylab("Latitude")
cruise_TA + dye_TA
Drifters
ggplot() +
geom_point(data = drifter_df, aes(lon, lat, shape = asset, color = log10(dye))) +
scale_color_viridis_c(option = "D",
name = "Dye (log ppb)",
trans = scales::log10_trans()) +
#geom_line() +
geom_point(size = 0.5) +
coord_quickmap() +
ggtitle("Drifter Dye") +
xlab("Longitude") +
ylab("Latitude")
Combined maps
Combined drifter, underway, and ctd (position only) data.
ggplot
Non-interactive
ggplot() +
geom_point(data = uw_dye_df,
aes(x = Longitude, y = Latitude, color = log10(Dye_ppb)),
size = .4) +
geom_point(data = drifter_df, aes(lon, lat, shape = asset, color = log10(dye))) +
scale_color_viridis_c(option = "D",
name = "Dye (log ppb)",
trans = scales::log10_trans()) +
geom_point() +
coord_quickmap()
Leaflet
Leaflet is slightly easier to work with, but much slower.
# Create a palette that maps factor levels to colors
pal_fac <- colorFactor(c("red", "navy"), domain = ctd_loc$in_out)
#pal_rho <- colorQuantile("magma", drifter_df[["Dye_ppb"]], n = 10)
pal_rho <- colorQuantile("magma", uw_dye_df[["sensor_ppb"]], n = 50)
leaflet(drifter_df) |>
addTiles() |>
addCircleMarkers(lng = ~lon,
lat = ~lat,
color = ~pal_rho(dye),
popup = ~asset,
radius = 2,
stroke = FALSE,
fillOpacity = 0.5) |>
addCircleMarkers(
data = uw_dye_df,
radius = 2,
stroke = FALSE,
fillOpacity = 0.5,
color = ~pal_rho(uw_dye_df[["Dye_ppb"]])) |>
addCircleMarkers(
data = ctd_loc,
color = ~pal_fac(in_out)
)Mapdeck map
Mapdeck uses WebGL for rendering, so is much faster with 70k points of underway data. Still figuring out color scales for long-tailed fluormeter data.
#m <- grDevices::colorRamp(c("blue", "white", "yellow"), bias = 6)(df[["Dye_ppb"]])
#pal <- col_quantile("viridis", df[["Dye_ppb"]], n = 30)(df[["Dye_ppb"]])
uw <- uw_dye_df |>
mutate(dye = log10(Dye_ppb))
dr <- drifter_df |>
mutate(dye = log10(dye))
key <- 'pk.eyJ1IjoiYmxvbmd3b3J0aCIsImEiOiJjbHI3NmtzMmYyZTBuMmluOGZ4dmhuaWFqIn0.eW7yUb_wqL7LQJOURH9IRQ'
mapdeck(token = key, style = mapdeck_style("light")) |>
add_scatterplot(
data = uw,
lat = "Latitude",
lon = "Longitude",
radius = 10,
fill_colour = "dye",
layer_id = "dye_layer",
palette = "viridis") |>
add_scatterplot(
data = dr,
lat = "lat",
lon = "lon",
radius = 10,
fill_colour = "dye",
layer_id = "drifter_layer",
palette = "viridis")